package fem2.observer;

import inf.jlinalg.ArrayVector;
import inf.jlinalg.IVector;
import math2.LinAlgUtilities;
import fem2.Debugger;
import fem2.Element;
import fem2.Model;
import fem2.element.StructuralElement;
import fem2.element.crack_sda.SDAElement;
import fem2.enu.StressType;

/**
 * observer to monitor maximum pricipal stress in the model (for SDAElement)
 * 
 * @author hbui
 * 
 */
public class StressObserverForSDAElement extends StressObserver {

	public StressObserverForSDAElement(Model m) {
		super(m);
	}

	@Override
	public void fireTimeStepFinished() {
		/*
		 * search for largest stress at Gauss point
		 */
		int ne = m.countElements();
		int dim = m.getMesh().getDimension();
		double o1Max = 0;
		int eIndex = 0;
		int gpIndex = 0;
		for (int i = 0; i < ne; i++) {
			StructuralElement e = (StructuralElement) m.getElement(i);
			if (e instanceof SDAElement) {
				int ng = e.getNumIntegrationPoints();
				for (int j = 0; j < ng; j++) {
					IVector sigma = e.getStress(j, StressType.CAUCHY_STRESS);

					// double[] pstress = LinAlgUtilities.spectralDecompose(dim,
					// sigma);
					// double o1 = pstress[0];

					IVector ps = new ArrayVector(dim);
					LinAlgUtilities.decomposePrincipalStress(dim, sigma, ps, null);
					double o1 = ps.get(dim - 1);

					if (o1Max == 0 || o1 > o1Max) {
						o1Max = o1;
						eIndex = i;
						gpIndex = j;
					}
				}
			}
		}
		Element emax = m.getElement(eIndex);
		// StructuralElement semax = (StructuralElement) emax;
		// IVector sigma = semax.getStress(gpIndex);
		Debugger.watch("Maximum sigma1 locate at " + emax.getElemType() + " SDAElement "
				+ emax.getMeshPart().getId() + " at Gauss point " + gpIndex + " : " + o1Max);
		// Debugger.watch("sigma = ", sigma);
		// for (int i = 0; i < emax.getNumGaussPoint(); i++) {
		// Debugger.watch(semax.getStress(i));
		// }
	}

	@Override
	public void fireFinished() {
	}

}
